%% Functions to plot intervals

function [fig, sortOrder] = plot_intervals(hat_vartheta, bar_tau, predict_se, highlight_indices, xtick, ytick)
    
    % Initialization
    N = length(bar_tau);

    % Calculate the recentered intervals
    re_lower_bounds = bar_tau - hat_vartheta - 1.282 * predict_se;
    re_upper_bounds = bar_tau - hat_vartheta + 1.282 * predict_se;
    
    % Sort intervals based on ascending lower_bounds, then by upper_bounds
    [~, sortOrder] = sortrows(re_lower_bounds);
    re_lower_bounds = re_lower_bounds(sortOrder);
    re_upper_bounds = re_upper_bounds(sortOrder);

    % Update highlight_indices based on sorted order
    highlight_indices = find(ismember(sortOrder, highlight_indices));

    fig = figure;
    set(fig, 'Units', 'inches', 'Position', [1 1 5.833 4.373]); 
    hold on;
    
    % Loop through each interval and plot bars, highlighting specified indices
    for i = 1:N
        % Check if the current index is in the highlight_indices list
        if ismember(i, highlight_indices)
            % Plot highlighted bar with a different color
            rectangle('Position', [i - 0.4, re_lower_bounds(i), 0.8, re_upper_bounds(i) - re_lower_bounds(i)], ...
                      'FaceColor', [0.6, 0.8, 1, 0.5], 'EdgeColor', 'none'); % Light blue
        else
            % Plot regular bar with pink semi-transparency
            rectangle('Position', [i - 0.4, re_lower_bounds(i), 0.8, re_upper_bounds(i) - re_lower_bounds(i)], ...
                      'FaceColor', [1, 0.6, 0.6, 0.5], 'EdgeColor', 'none'); % Pink
        end
    end
    
    yline(0, 'k--', 'LineWidth', 1);
    yticks(ytick);
    ylim([ytick(1), ytick(end)]);
    xticks(xtick);
    xlim([xtick(1), xtick(end)]);
    set(gca, 'FontSize', 18);
   
    hold off;

end
